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' Abstract. We present the stability analysis of the dark states in the adiabatic 



passage for the linear and non-linear lambda and tripod systems in the presence of 
amplitude damping (losses). We perform an analytic evaluation of the real parts of 
eigenvalues of the corresponding Jacobians, the non-zero eigenvalues of which are found 
from the quadratic characteristic equations, as well as by the corresponding numerical 
simulations. For non-linear systems, we evaluate the Jacobians at the dark states. 
Similarly to the linear systems, here we also find the non-zero eigenvalues from the 
characteristic quadratic equations. We reveal a common property of all the considered 
systems showing that the evolution of the real parts of eigenvalues can be split into 
three stages. In each of them the evolution of the stimulated Raman adiabatic passage 
(STIRAP) is characterized by different effective dimension. This results in a possible 
adiabatic reduction of one or two degrees of freedom. 
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1. Introduction 

Over the last couple of decades there has been a continuing interest in the stimulated 
Raman adiabatic passage (STIRAP) [1-5]. The simplest situation is the adiabatic 
passage in a linear lambda system [6-12] containing a single dark (uncoupled) state 
which is immune to the atom-light coupling. If atomic initial and final states are the 
ground states representing the dark states of the system, the atom can be transferred 
between these two states by slowly changing the relative intensity of the laser pulses. 
When the adiabatic passage is slow enough, the excited state is only slightly populated 
and thus the losses are minimum. The analysis has been extended for the STIRAP 
process in the tripod system characterized by two dark states [13-15]. This enables 
to create a quantum superposition of metastable states out of a single initial state 
in a robust and coherent way [16, 17]. The schemes involving more atomic and 
molecular levels were also proposed for creation a superposition of states [18] as well 
as for experimental control of excitation flow [19]. Recently the treatment was further 
extended to the non-linear lambda [20-31] and tripod [32] schemes. 

Usually the STIRAP is based on the adiabatic approximation. However, one has 
to distinguish between the adiabatic approximation and the adiabatic reduction of 
dynamic systems. In quantum mechanics, a closed quantum system is said to undergo 
adiabatic dynamics if its Hilbert space can be decomposed into decoupled Schrodinger 
eigenspaces with distinct, time-continuous, and non-crossing instantaneous eigenvalues 
of Hamiltonian. On the other hand an open quantum system is said to undergo adiabatic 
dynamics if its Hilbert-Schmidt space can be decomposed into decoupled Lindblad- 
Jordan eigenspaces with distinct, time-continuous, and non-crossing instantaneous 
eigenvalues of the Lindblad superoperator [33]. The system is called adiabatically 
approximated if the error term in the Schrodinger equation (for closed systems) or in the 
master equation (for open systems) is much less than the diagonalysed part; i.e. one may 
neglect the non-diagonal terms in order to get an adiabatically approximated version 
of the system. Note that in the presence of fast driven oscillations some additional 
conditions (in addition to the slowness of the evolution of the Hamiltonian) have to be 
imposed [34]. 

Another procedure is the adiabatic elimination of decaying degrees of freedom. 
It is related to the dynamic systems in which some degrees of freedom may decay. In 
accordance with this definition, these degrees of freedom may be adiabatically eliminated 
by solving the corresponding algebraic equations; the r.h.s. of decaying equations are set 
to be equal to zero. Consequently, one obtains the dependencies of decaying variables 
on the remaining ones. If there are several decaying variables, one may eliminate them 
one by one, starting from the fastest variable, and finishing with the slowest one. Such 
a procedure can be found e.g. in the book by Haken [35]. 

The aim of the present work is to perform a stability analysis of the dark states 
in the adiabatic passage for the linear and non-linear lambda and tripod systems. The 
analysis sets limits to the adiabatic reduction of the systems. Moreover, we have revealed 
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that in all the considered systems, the stability properties of the dark states are similar, 
namely, there are three time intervals with different number of the negative real parts of 
the Jacobians. This suggests that the corresponding linear and non-linear systems have 
equal possibilities of adiabatic reduction. Although the linear lambda [1-12] and tripod 
systems [13-17, 36] have been substantially studied in the literature, here we apply our 
treatment also to these linear setups in order to facilitate the subsequent analysis of the 
non-linear systems. 

The non-linear lambda system can be realized in the Bose-Einstein condensates 
(BEC) via photoassociation (PA) from a dissociated (quasicontinuum) atomic state 
to the ground molecular state in the presence of the intermediate molecular state 
[20,21]. The aim of the quantum control is to transfer the whole population from the 
dissociated atomic to the ground molecular state. One thus creates ultracold molecules 
by associating cold atoms [29,37]. In this case the dark state is a generalisation of that 
for the linear systems. In the linear case, the dark state is defined as a superposition 
of the initial and target ground atomic states which corresponds to the zero eigenergy 
of the system Hamiltonian. The same dark state may be also defined as a steady state 
solution of the Schrodinger equation. If we consider the nonlinear system, we can again 
define the dark state as a steady state [29] of the equations of motion following from 
the Heisenberg equation. In the nonlinear systems, the behaviour of the dark state 
reproduces that in the linear systems although the superposition is missing now. We 
thus get the nonlinear version of STIRAP: the entire population is distributed among 
the steady state probability amplitudes of the initial atomic and the target molecular 
state. At the beginning the whole population is atomic, whereas at the end it is in the 
ground molecular state. 

Differently from the traditional STIRAP in an atomic lambda system, the atom- 
molecule STIRAP contains nonlinearities originating from the conversion process of 
atoms to molecules, as well as from the interparticle interactions described by the non- 
linear mean-field contributions. The existence of such nonlinearities makes it difficult 
to analyse the adiabaticity of the atom-molecule conversion systems due to the absence 
of the superposition principle. In the STIRAP, the linear instability could make the 
quantum evolution deviate from the dark state rapidly even in adiabatic limit [22]. 
Therefore, it is important to avoid such an instability for the efficiency of the STIRAP. 

The non-linear tripod system can be realized in the PA with two target states 
involved. Specifically, one may consider the atom-molecule transition in ultracold 
quantum gases via PA. It was first considered in [32] where the second order dynamic 
system was derived that parametrizes the solution evolving on the dark state manifold. 
However the stability of the solution moving along the manifold was not considered. 
Therefore we shall check the stability of this solution, i.e. see whether the nearby 
solutions are attracted back to this manifold. 

The adiabatic theory for non-linear quantum systems was first discussed by Liu et 
al. [38] who obtained the adiabatic conditions and adiabatic invariants by representing 
the non-linear Heisenberg equation in terms of an effective classical Hamiltonian. Pu 
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et al. [23] and Ling et al. [24] extended such an adiabatic theory to the atom-dimer 
conversion system by linking the nonadiabaticity with the population growth in the 
collective excitations of the dark state. Specifically it was shown that a passage 
is adiabatic if the solution remains in a close proximity to the dark state. Itin 
and Watanabe [25] presented an improved adiabatic condition by applying methods 
of the classical Hamiltonian dynamics. The atom-molecule dark-state technique in 
the STIRAP was theoretically generalized to create more complex homonuclear or 
heteronuclear molecule trimers or tetramers [39-42]. 

An important issue is the instability and the adiabatic property of the dark state 
in such complex systems. For example, the dynamics of a non-linear lambda system 
describing BEC of atoms and diatomic molecules was studied and a model of the dark 
state with collisional interactions was investigated [26, 27]. It was shown that non-linear 
instabilities can be used for precise determination of the scattering lengths. On the 
other hand, the transfer of atoms to molecules via STIRAP is robust with respect to 
detunings, x 3 nonlinearities, and small asymmetries between the peak strengths of the 
two Raman lasers [27]. The complete conversion is destroyed by spatial effects unless 
the time scale of the coupling is much faster than the pulse duration. In addition, a 
set of robust and efficient techniques has been introduced [43] to coherently manipulate 
and transport neutral atoms based on three-level atom optics (TLAO). 

It is to be noted that the dynamics of an adiabatic sweep through a Feshbach 
resonance was studied [44] in a quantum gas of fermionic atoms. An interesting 
application of BEC is an atom diode with a directed motion of atoms [45]. Another 
example of BEC was presented in [46] where it was shown that the two-colour PA of 
fermionic atoms into bosonic molecules via a dark-state transition results in a significant 
reduction of the group velocity of the photoassociation field. This is similar to the 
Electromagnetically induced transparency (EIT) in atomic systems characterized by 
the three-levels of the lambda type. In addition the coupled non-linear Schrodinger 
equations have been considered [47] to describe the atomic BECs interacting with the 
molecular condensates through the STIRAP loaded in an external potential. The results 
have shown that there is a class of external potentials where the exact dark solutions 
can be formed. In [48] it was shown that it is possible to perform qubit rotations by 
STIRAP, and proposed a rotation procedure in which the resulting state corresponds 
to a rotation of the qubit, with the axis and angle of rotation determined uniquely by 
the parameters of the laser fields. 

A relevant tool for studying the adiabaticity is the adiabatic fidelity. It indicates 
how far is the current solution of the system from the dark state. Meng et al. have 
generalized the definition of fidelity for the non-linear system [28]. They have studied 
the dynamics and adiabaticity of the population transfer for atom-molecule three-level 
lambda system on a STIRAP. It was also discussed how to achieve higher adiabatic 
fidelity for the dark state through optimizing the external parameters of the STIRAP. 
In the subsequent paper [49] Meng et al. have used the same definition of adiabatic 
fidelity in order to discuss the adiabatic evolution of the dark state in a non-linear atom- 
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trimer conversion due to a STIRAP. It is to be noted that Ivanov and Vitanov have 
recently proposed novel high-fidelity composite pulses for the most important single- 
qubit operations [50]. 

In this work, we analyse the problem of reducing the dimension (simplifying) in the 
linear/non-linear three- and four-level models. This procedure is called the adiabatic 
reduction, and its validity is closely related with the theory of adiabaticity discussed 
above. The exact three- or four-level system may be adiabatically reduced to a system 
with lower dimension. The question that arises is how many dynamic variables can be 
eliminated? In other words, what is the effective dimension of the system? The answer 
lies in the eigenvalues of the Jacobian computed at the dark state. The zero real parts 
of eigenvalues mean that in some directions the nearby solutions are behaving neutrally 
in respect to the dark state. The negative real parts in turn mean that some directions 
are stable, and the nearby solutions converge towards the dark state. Therefore we 
conclude that the number of negative real parts dictates the number of variables that 
can be adiabatically eliminated (see e.g. [35]). On the other hand, the number of zero 
real parts yields the effective dimension of the system. Note that we find the non-zero 
eigenvalues analytically from quadratic characteristic equations. 

One of the central issues in our work is the presence of dissipation in all the 
considered systems. The non-zero losses make the adiabatic reduction easier to 
implement since the term of losses acts as a "controller" that attracts the nearby 
solutions towards the dark state. However, Vitanov and Stenholm have demonstrated 
that the losses cause also the decrease of transfer efficiency to the target state [8]. 
This decrease can be circumvented by higher pulse areas since the range of decay rates 
over which the transfer efficiency remains high, has been found to be proportional to 
the squared pulse area (see (10) in [8]). In the subsequent developments the effect of 
spontaneous emission on the population transfer efficiency in STIRAP was explored for 
the linear lambda [10, 11] and tripod [15] systems. In addition, Renzoni et al. [51] have 
considered the coherent-population-trapping (CPT) phenomenon in a thermal sodium 
atomic beam. It was demonstrated that CPT may be realized on those open transitions 
with an efficiency decreasing with the amount of spontaneous emission towards external 
states. On the other hand, here we concentrate on the stability issues of the linear and 
non-linear lambda and tripod systems with the losses. 

The paper is organized as follows. In next two sections we shall consider the stability 
of the linear lambda and tripod systems with losses. In sections 4 and 5 the analysis is 
extended to the non-linear lambda and tripod systems. In section 6 we discuss the role 
of the one-photon detuning followed by the conclusions in section 7. 

2. The linear lambda system 

In this section we shall provide a summary on the STIRAP in the linear lambda system 
with losses studied in [8,10,11] followed by the stability analysis of the system. The 
three- level lambda system is shown in figure 1. The excited state |e) is coupled to two 
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ground states |a) and |g) with the coupling strengths denoted as fi p and fid to form a 
lambda scheme. The Hamiltonian of such system reads: 

H = -h(A + i T )| e ><e| + |[O p |a)(e| + O d |g><e| + H.c.]. (1) 




Figure 1. Three-level system coupled by two lasers. f2 p and fid are the Rabi 
frequencies for the pump and damp laser, A is the one photon detuning, and 7 is 
the loss rate. 

Note that in this Hamiltonian besides the real-valued one-photon detuning A there 
is the imaginary term 7 representing the losses. Denoting the amplitudes of the three- 
level state as ip ei ipa., and ipg respectively, we get the Schrodinger equation: 

k/> a = tt p 1p e , (2) 
k/>e = - (A + i7)^ e + fi p ^ a + fi d ?/>g, (3) 

k/> g = n d ip e . (4) 
The normalization reads: 

l^a(t)| 2 + |^(t)| 2 + |^e(t)| 2 <l, (5) 

where equality holds for initial time. Because of losses (7 > 0) the total normalization 
will be slightly reduced (for t > 0) during the transfer of population through the excited 
level. This property of the total population is assumed throughout the paper. 
We take the laser pulses to be Gaussian 

(t-td) 2 



fid = f2 exp 



fl exp 



J^2 
J^2 



21 



(6) 
(7) 



Here the pulses are centered at t d and t p , respectively, T is their width, and Q is their 
amplitude. 

The third order system (2,3,4) may be rewritten in a matrix form 

* = -iH^f = (8) 
where ^ = [if) a , i[) e , ip g ] T is the vector of the state of the system, and 

/ fi p \ 
// fi p -(A + i 7 ) o, (9) 
fi d 



V 







/ 
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is the corresponding Hamiltonian. The matrix A = —iH is the Jacobian of the system. 
If the Hamiltonian possesses eigenvalues u, then the eigenvalues of the Jacobian are 
defined by A = — iui. Note that the real parts of A determine the stability of the fixed 
point at the origin. 

We now find these eigenvalues, more specifically their real parts. The eigenvalues 
of Hamiltonian satisfy the characteristic equation 

det \\H-Iu\\ = 0, (10) 

with / denoting the unit matrix. Expanding the corresponding third-order determinant 
one finds that one eigenvalue is always zero: 

wi = 0. (11) 

The other two eigenvalues satisfy the quadratic equation: 

u 2 + (A + i 7 v - (n 2 + n 2 p ) = o. (12) 

Thus the two eigenvalues satisfy 

u 2 + u 3 = - [A + i 7 ], (13) 
U2U3 = -[Q 2 d + Q 2 p ]. (14) 

Assuming A = 0, the solutions of (12) read: 

W2 ,3 = {-i 7 ± [-7 2 + 4(n2 + n 2 p )] 1/2 }/2. (is) 

For t — > ±oo the Rabi frequencies go to zero, and it follows from (13,14) that 
u>2 = 0,w 3 = — i 7 . Calling on (11), we can write 

Ai, 2 = 0, A 3 = -7, (16) 

or Re(Ai j2 ) = 0, and Re(A 3 ) = —7 for t — > ±00. 

For finite times there is a region where the Rabi frequencies are large enough, so 
that the discriminant is positive in (15): D = — 7 2 + 4(fi^ + Q^) > 0. Such a situation 
occurs in a certain interval t± < t < £2, and from (15) we get 

u; 2 ,3= -i 7 /2±v^/2, (17) 
A 2 , 3 = -7/2Tiv^D/2. (18) 

The first eigenvalue is Ai = u\ — 0. The boundaries ti,t 2 are solutions of D(i) = in 
respect to time. 

Hence, in the interval t± <t < t 2 , the real parts are Re(Ai) = 0, Re^^) = —7/2. 
We may adiabatically reduce the dimension of this system but we first transform 
its variables. We change the bare variables i[) a , ip g to the bright ips and dark ipn one: 

^ B = (^a + n d ^ g )M (19) 

= (fid^a " (20) 

where 

^ = (^ + ^d) 1/2 - (21) 
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Denoting £ p = Q p /Q, £d = ^d/^, and performing some operations, we obtain the 
following equations for the new variables: 

k/> B = onfo + TQij e , (22) 

i^e = -T(A + i 7 )^c + ™/>B, (23) 

iip D = a*ifo, (24) 

where we have made the time dimensionless by substituting t/T — > t. Here a = 
i(£p£d _ £d£ P ) is a dimensionless parameter in which derivatives £ p , £d are taken in respect 
to dimensionless time. 

We now adiabatically eliminate the amplitude ip e by setting ip e = 0. From (23) we 

get 



(25) 



A + 17 

Inserting this result in (22,24) we obtain 

TVL 2 

iip B = aip D + ifo, (26) 

A + 17 

iipD = a*ip B - (27) 

We solve this system to find the dynamics of ipB, ipD, and from (25) we find the dynamics 
of ip e . 

The system (26,27) may be also adiabatically reduced. We now set ipB = 0, and 
solve (26) for ip B : 

a(A + i 7 ) 



TO? 



-ifo- 



Inserting this result in (27), we get a first order dynamic system 

|a| 2 (A + i 7 ) 



i^D 



TQ 2 



The Hamiltonian for the 2D reduced system (26,27) reads: 

TO 2 



if. 



2D 



, . a 
A+17 

a* 



(28) 



(29) 



(30) 



Solving the characteristic equation for this Hamiltonian, and assuming A = 0, one 
obtains two eigenvalues: 



W12 



TO 2 

7 
Tfl 2 



± 



THl 4 
— 

T 2 n 4 
— 



± 



7 



- A\a\ 



- A\a\ 



1/2 



1/2- 



(31) 



(32) 



Here Ai j2 are the eigenvalues of the corresponding Jacobian defined as A 2D = —iH 2 b. 
There are two regimes of evolution of the eigenvalues corresponding to D 2 < and 
D 2 > 0, where D 2 = - A\a\ 2 is the discriminant in (31,32). 



Stability of lambda and tripod systems 



9 



Hence, we have the eigenvalues for all the three versions of the linear lambda system. 
For the exact 3D system, they are given by Ai = U\ = and (17,18). For the reduced 2D 
system, they read (31,32). Finally, a trivial single eigenvalue for the ID system follows 
from (29): 

|o;| 2 7 

Ai = -i^i = - w . (33) 

In figure 2 we show the dynamics of the real parts of eigenvalues for all three cases. 

As the solid lines in figure 2(a) show, for t < t 1 and t > t 2 there are two different 
non-zero branches for the 3D case: the upper branch determines the slow decay of 
the bright state, whereas the lower branch causes the fast decay of the excited one. 
Therefore, we may eliminate the excited state, and cannot do this with the bright one. 
In the middle of the process, where t\ < t < t 2 , both branches become degenerate with 
real part of eigenvalues equal to —7/2. We thus may reduce the both states, excited 
and bright. During the whole evolution, one real part remains exactly zero showing that 
the dark state does not experience any losses. This fact indicates that the process is 
adiabatic. 

The short dash dotted and short dash lines show the dynamics of two real parts 
for the 2D reduced system. Although the time moments [(ti,t 2 ) = (1.24,5.56)] at 
which the discriminant D 2 is zero differ from the times [(^1,^2) = (1.49,5.32)], where 
the discriminant D goes to zero, these differences are small compared to the interval 
(t 2 — ti). In the time intervals of t < t\ and t > t 2 the both branches are close to zero 
indicating that neither bright nor dark state may be eliminated. The process is here 
again 2D. In contrast, in the range of time, where t\ < t < t 2 , the decay rate of the 
bright state is very large compared to that of the dark state. It is well seen in the figure 
2(b). Therefore the process is again ID as for the exact case. 

And lastly, the dotted line shows the dynamics of the real part for the ID reduced 
system. One may clearly see that for t < t± and t > t 2 the decay of the dark state is 
rapid compared to that for t\ < t < t 2 . The fastness of decay for t < t\ and t > t 2 
is seen expressively in the figure 2(b). This leads to conclusion that the ID reduced 
system is appropriate only in the time interval of ti < t <t 2 . 

Exactly as for the 3D system, the evolution of the system in the 2D (ID) 
approximations is adiabatic for the whole time of integration (in the time interval 
ti < t < t 2 ) since for both cases the decay rate of the dark state is very small compared 
with that of the bright state. 

3. The linear tripod 

The STIRAP process in the linear tripod scheme without dissipation was first considered 
by Unanyan et ai [13, 14]. Here we outline this scheme in which the dissipation is also 
included. Afterwards we perform the linear analysis of this system. 

Consider the four-level system schematically shown in figure 3. The excited state 
|e) is coupled to three ground states |a), |gl), and |g2) with the coupling strengths 
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Figure 2. (a) Dynamics of real parts of eigenvalues of the Jacobian for linear lambda 
system, computed for all three cases (3D, 2D, and ID) by (17,18), (31,32) and (33), 
respectively. The dashed vertical lines set the boundaries for the ID and 2D processes. 
Here t\ = 1.49 and ti = 5.32. The solid lines are the real parts for the 3D, the short 
dash dotted and short dash line correspond to the 2D case, and the dotted line is the 
real part for the ID reduced system; (b) the dynamics of real parts for the 2D and ID 
cases in extended vertical scale. The lines are chosen in the same way as in (a). The 
parameters are as follows: A = 0, 7 = 2.0, fio = 10.0, t p = 3.8, id = 3.0, and T = 1.0. 



denoted as Q p , f2di, and f2<j2, respectively. Here p stands for pump and d stands for the 
damp. The four-level Hamiltonian reads: 



H = -h{A + ij)\e)(e\ + f [fi p |a)(e| + fi dl |gl)(e| + 
+n d2 |g2)(e|+H.c.]. 



(34) 




Figure 3. Four-level system coupled by three lasers. f2 p and fidi, ^d2 are the Rabi 
frequencies for the pump, and damp lasers, A is the one photon detuning, and 7 is the 
loss rate. 

Denoting the amplitudes as i/j e , ^> a , if) s i, and ipg2, the Schrodinger equation reads: 
iip a = Qpip e , (35) 

i^ e = - (A + i7)^ e + fip^a + ^dl^gl + 

+ ^d2^g2, (36) 
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|^a(t)| 2 + |^ g i(t)| 2 + \ri>#(t)\ 2 + \Mt)\ 2 < 1, 

where equality holds for initial time. 
The Rabi frequencies are given by 

(t-t P r 



11 

(37) 
(38) 

(39) 



f2 p = f2 exp 
Q d i = Kifl exp 
il d2 = -^2^0 exp 



Ji2 

(t - tdl) 2 
(t - t d2 ) 2 

J"2 



(40) 
(41) 
(42) 



H 



V 



(43) 



Here the pulses are centered at t p , t dl and t d2 , respectively. T is the width of the pulses, 
Ki, K 2 determine the amplitudes for the damp pulses, and Q is the amplitude for the 
pump pulse. 

The system (35,36,37,38) can be written in the form of (8) with the state vector 
^ = ip e , ipgi, i J g2] T , and Hamiltonian 

/ tt p \ 

n p -(A + i 7 ) n dl n d2 

On 

o n d2 o o j 

The matrix A = —\H is again the Jacobian of the system. Here the relation A = — ioo 
holds. Solving the eigenvalues problem (10) for the linear tripod, we obtain two zero 
eigenvalues, 

wi,2 = 0. (44) 
The other two eigenvalues can be found from quadratic equation 

oo 2 + (A + i 7 ) w - (n 2 dl + n 2 d2 + n 2 p ) = o. (45) 

The eigenvalues cu 3i 4 must satisfy 

+ = - (A + i7), (46) 

u 3 u 4 = -(n 2 dl + n 2 d2 + n 2 p ). (47) 

We again assume that A = 0, thus obtaining the following solutions: 

c 3 ,4 = {-i7 ± [-7 2 + 4^ + Q 2 d2 + n 2 p )] 1/2 }/2. (48) 

(See also (7) in [14]). The Rabi frequencies are chosen in the form of Gaussian pulses. 
For t — > ±oo the Rabi frequencies tend to zero, and from (46,47) it follows that cu 3 = 0, 
004 = —vy. Hence, for t — > ±00 we have 

Ai,2,3 = 0, A 4 = -7- (49) 
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However between these two infinite times there is a region where the Rabi frequencies are 
large enough, and the discriminant in (48) is positive [D = —7 2 + 4(1^ + 1^2 + 1^) > 0]. 
Such a situation takes place in the interval t\ <t<t2, and from (48) we get 

w 3 ,4 = - i-y/2 ± V3/2, (50) 
A 3 ,4 = - 7/2 T iVD/2. (51) 

The first two eigenvalues are Ai j2 = wi,2 = 0. The boundaries ti,t 2 are defined in the 
same way as in section 3. 

Hence, in the interval t± < t < t 2 , the real parts are Re(Ai j 2) = 0, Re^^) = —7/2. 

Exactly as in the previous section, we first transform the variables from the bare 
states to one bright and two dark states: 

(V> a , A, lf>gl,M ->■ (V'B, ^e, ifol, ^2) ■ (52) 

Parametrizing the Rabi frequencies as 

f2 d2 = Q sin (0), (53) 

Q p = ficos(0)sin(0), (54) 

Qdi = & cos (0) cos (0), (55) 
we may write down the amplitudes of one bright and two dark states: 

ipB = cos (0) sin (0)-0 a + cos (0) cos (G)V'gl + 

+ sin(0)0> g2 , (56) 

ip m = cos (0)^ a - sin (0)^ g i, (57) 

■0D2 = sin (0) sin (0)-0 a + sin (0) cos (0)^ g i - 

- cos (0)^g2- (58) 

After renormalizing the time (t/T — > t) and some rearrangements, we derive the 
following dynamic system for these variables: 



dt 



\ ^D2 J 






TQ 


«13 


a 14 


Ttt 


— T(A + i 7 ) 








* 

"13 








«34 


"14 





"34 






\ 



/ "0B \ 
\ ^D2 J 



(59) 



Here 



"13 = i(£n6tt +6363), 

"14 = 1(^11^41 + ^13^43 + C14C44), 
«34 = 1(^31^41 + ^33^43), 

and the coefficients £y are defined by the matrix 

^ cos (0) sin (0) cos (0) cos (0) 

1 

cos(0) -sin (6) 

y sin (0) sin (0) sin (0) cos (0) 



(60) 
(61) 
(62) 



lie 



ij 1 



sin (0) \ 



-COS(0) J 



(63) 
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Note that this matrix realizes the transformation (52) (see also (56,57,58)). 

We now reduce the dimension of system (59) in two steps. For the first step, we 
eliminate the excited state by setting 

4 = 0. (64) 
Solving the equation 

fi^B - (A + i 7 )^ e = 0, (65) 
(see the second row in (59)) for ip e , we get 



(66) 



A + i 7 

Inserting this result in (59), we obtain the following three-dimensional dynamic system: 



dt 



( 






/ 




ifol 






\ 


^D2 


J 


\ 



A+i7 



«13 "14 



Ot 



13 



\ 


/ 


4>B 


\ 










/ 


V 


^D2 


/ 



(67) 



a 34 

• : 1 ^34 

We can solve this system to find the evolution of ipB, V'di, i*D2, and to determine the 
dynamics of the excited state using (66). 

For the second step, we eliminate the bright state, i.e. we set 

= 0, (68) 

in (67). After solving the equation 

TQ 2 



-IpB + ai3^Dl + «14^D2 = 



A + i 7 

for ipB (see the first row in (67)), we find 

A + i 7 , 



TQ 2 



■(«13^D1 + «14^D2)- 



Inserting this result in (67), we obtain the following second order system: 





The matrix of this system may be split into two parts: 
H 2 



2 TQ 2 2 



where 



and 



(o) 



H. 



n 2 — 



a 34 

«34 



I«i3| a 13 au 

ttl4«13 |«14| 2 



(69) 



(70) 



(71) 



(72) 



(73) 



(74) 
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When Q 2 is relatively large, one can neglect the influence of H^p and write approximately 
H 2 — H^. Thus one arrives at a simple second order system 








a 



34 




(75) 



which is equivalent to (27) of [13]. However, this system is not an adiabatically 
reduced version of the system (35,36,37,38). Actually, it determines the solution moving 
exactly on the dark state manifold comprising the two degenerate states \D1) and \D2) 
(this statement can be confirmed by applying the approach of [32]). Comparing (72) 
with the corresponding result for the linear lambda system, (29), we note that the 
both Hamiltonians contain the characteristic time scale T (the pulse width) in the 
denominators. In the case of (72) this time scale is involved only with the correction 
Hamiltonian H%~\ and it is absent in since it corresponds to a zero on the r.h.s. of 
(29). The correction in (72) thus corresponds to the r.h.s. of (29). 

Figure 4 shows the dynamics of the real parts of the eigenvalues of the Jacobians 
for various dimensions. In figure 4(a) the solid lines correspond to the exact 4D case. 
For t < ti and t > t 2 the lower branch, for which ReA ~ —7, causes the decay of the 
excited state. The upper branch describing the decay rate for the bright state is small 
compared to the decay rate of the excited state in these time intervals indicating that 
the excited state may be adiabatically eliminated, whereas the bright state must be left. 
In the middle of the passage (ti < t < t 2 ) the decay rates of the excited and bright 
states become degenerate and equal to —7/2, thus enabling to eliminate them both. 
The eigenvalues for the dark states remain both zeros for the whole evolution of the 
system demonstrating that the process is adiabatic since the degenerate dark state does 
not loose its population. 

The dash dotted line in figure 4 displays the dynamics of the decay rate for the 
bright state in the 3D system. In figure 4(a) it falls down (grows up) just after t\ (just 
before t 2 ). In figure 4(b) the same dynamics is shown in extended vertical scale. From 
both figures, 4(a) and 4(b), one may conclude that in the 3D reduced system the bright 
state may be eliminated for t± < t < t 2 , and it should be preserved for t < t\ and 
t > t 2 , since its decay rate in the latter case is much less than in the former one. The 
two eigenvalues corresponding to degenerate dark states remain exactly (or almost) zero 
for the 3D system. They are small compared with the decay rate of the bright state. 
Therefore these states may not be eliminated for the whole time of evolution. 

The dotted and short dotted lines in figure 4 show the dynamics of the real parts 
for the 2D reduced system. One of them grows rapidly before t\ and converges to zero, 
whereas the other one is first zero and then decays rapidly after t 2 . Such a situation 
suggests that the range of applicability of the 2D system should be wider than the time 
interval t± < t < t 2 since the two degenerate dark states survive (do not decay) when the 
decay rate is relatively small. However, such a conclusion would be true if we distributed 
the whole population among the degenerate states at the end of the rapid growth of 
the first real part. But if the bright state had some initial population, it could not be 
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neglected since its real part is close to zero before t\ and after t 2 . 

We should also note here that for the 2D reduced system, the decay rate for the 
second dark state grows rapidly from negative values to zero (the decay rate for the first 
dark state decreases rapidly from zero to negative values) in the time intervals t < 2.0 
(t > 18.0), i.e. outside of the figure 4. But we do not need to take these events into 
account since we are interested in the time interval where both decay rates are close 
to zero. This interval is determined by the growing (decreasing) decay rates which are 
plotted in the figure 4. 

Exactly as in the previous section, we here can also conclude that the adiabaticity 
is preserved also for reduced systems, since for the 3D (2D) approximations the dark 
states manifold does not loose its population for the whole time of evolution (in the 
time interval that is wider than t\ < t < t 2 ). 




Figure 4. (a) Dynamics of real parts of eigenvalues of the Jacobian for linear tripod, 
computed for all three cases (4D, 3D and 2D) by [(48), (44)], from (67) and (72), 
respectively. The dashed vertical lines set the boundaries for the 2D and 3D processes. 
Here t\ = 6.12 and i 2 = 12.73. The solid lines are the real parts for the 4D, the short 
dash dotted and short dash line correspond to 3D case, and the dotted and the short 
dotted lines are the real parts for the 2D reduced system; (b) the dynamics of real parts 
for the 3D and 2D cases in extended vertical scale. The lines are chosen in the same 
way as in (a). The parameters are as follows: A = 0, 7 = 2.0, Qq — 60.0, t p = 10.7, 
i d i = 10.0, t d2 = 8.5, K x = 0.75, K 2 = 5.0, and T = 1.0. 

Figure 5 displays the results of numerical computations. In figure 5(a) we have 
presented the dynamics of populations of the bright (dashed and dotted lines), and 
excited (dashed-dotted and short-dotted lines) states. The dashed and dotted lines 
correspond to solutions of exact equations (35,36,37,38). The dashed-dotted and short- 
dotted lines are obtained from (67) (the bright state), and from (66) (the excited state). 
One can see that the populations of these states remain relatively small (P c , Pb — 10~ 3 ). 
In figure 5(b) the dynamics of coupling strengths is shown. Figures 5(c,d) display 
the dynamics of populations of the degenerate dark states. Figure 5(c) presents the 
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solutions of exact (35,36,37,38) (thin solis line) and that of adiabatically reduced system 
(67) (thick solid line). The exact and approximate solutions are in good quantitative 
agreement (we do not distinguish them in the present graph). Both systems were 
integrated in the whole time range (t G [0.0,20.0]). In figure 5(d) we compare the 
dynamics for exact (35,36,37,38) (thin solid line) with those of adiabatically simplified 
(71) (thick solid line). 




5 10 15 20 ,„ 5 10 15 20 
(C) t (d) t 



Figure 5. (a) The dynamics of populations 10 3 P e ,10 3 PB that are computed from 
(35,36,37,38) (dashed/dash-dotted lines), and from (67), (66) (dottcd/short-dotted 
lines); (b) the dynamics of Gaussian pulses computed by (40,41,42); (c,d) the dynamics 
of populations of dark states computed in (c) by (35,36,37,38) (thin solid line), and 
by (67) (thick solid line); in (d) the thin solid line represents again the exact solution 
of (35,36,37,38), and thick solid line is obtained using the approximation from (71) 
computed in the range of t\ <t<t 2 - The parameters are the same as in figure 4. The 
initial conditions are Va(0) = 1, V'gi(O) = ^2(0) = V'e(O) = 0. 



In figure 6 we compare the exactness of various approximations for the non-linear 
tripod. We see that the system (75) is in a good quantitative agreement with the exact 
system, but it is worse than (67) in the range t > t±. Whereas the 2D approximation 
(71) coincides with the solution of (67) in the interval t\ < t < t 2 almost identically. At 
the end of this interval, just before t 2 , the solution of (71) slightly deviates from that 
of (67). This is due to the fact that the magnitude of element ^2 1 a 13 | 2 contributed by 
the matrix in (72) becomes large compared with the magnitude of the element 0:34 
of the matrix Hf \ (73). The remaining elements of H^p are much less than | 1 2 just 
before £2- 

We should also note that in figure 6 the solution of "non-exact" 2D system (73) 
coincides almost identically with the exact 4D solution in the range of t < t\. However, 
such a situation takes place when the initial conditions are very close to the dark state. 
If we pushed them away from the dark state, the result of (73) would become worse 
than that of the 3D system (67) (not shown here). The system thus remains effectively 
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3D in the range of t < t\. 

Similarly it can be numerically verified that the solution of (71) is not worse than 
that of (67) just before t<i if one pushes the initial conditions further away from the dark 
state manifold. Therefore the system remains 2D in the whole interval t\ < t < t 2 . 



0.0000 



S -0.0005 



ra □ 

Q_ 



a. -0.0010 




Figure 6. The differences between the exact and approximated degenerate dark states; 
the short dotted and dotted lines show the differences of Pdi — Pdi an( ^ ft>2 — -P§2 f° r 
the 3D system (67), respectively; the dash dotted and short dash dotted lines display 
the corresponding differences for Pdi and Pd2 for the 2D system (75), respectively; 
the solid and dashed lines show the differences for Pdi and Pd2 for the adiabatically 
reduced system (71), respectively (integrated only in the range of t\ < t < £2)- The 
vertical dashed lines set the boundaries for the 3D and 2D processes. The parameters 
and initial conditions are the same as in figure 5. 



4. The non-linear lambda system 

The three-level non-linear Hamiltonian for the non-linear lambda system (see the figure 
1) reads: 

+H.c). 1 ' 

Here ^ a ,^ (a = a, e,g) are the bosonic annihilation and creation operators for state 
a, respectively. When the number of particles is much larger than the unity, the boson 
operators are replaced by c numbers tp a ,ipa (the meanfield treatment [29,37]) which 
obey the following Heisenberg equations 

k/> a = fip^e, (77) 

i4 = - (A + i 7 )^ c + hl^l + ifi d ^ g , (78) 
i4 = ^ d Vv (79) 
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The normalization reads: 

|^ a (t)| 2 + 2[|^(t)| 2 + |^(t)| 2 ]<l, (80) 

where the equality holds for the initial time. 

The nonlinearity enters here in the coupling induced by the pump field: it couples a 
particle in the state e with a pair of particles in the state a. Such a nonlinear couplig is 
encountered in second harmonic generation in nonlinear optics (where a represents the 
fundamental photon and e its second harmonics), as well as in the photoassociation of 
atoms into diatomic molecules [20,21,29,37], where a represents an atomic state, while 
e and g are excited and ground diatomic molecular states, respectively 

As for the linear lambda system, we take the Gaussian pulses given by (6,7). 

Similarly as in two previous sections, we here define the state vector * = 
[■0a, V'ej i^gV- The dynamic system (77,78,79) can be rewritten in the vector form: 

4* = /W, (81) 

where / is the vector of (generally) non-linear functions on the r.h.s. of system 
(77,78,79). The steady state solution ^o(t) of this system represents the dark state 
which is obtained by solving 

/(*o) = 0. (82) 

For the non-linear lambda system (77,78,79) the dark state reads ([29,37]): 

2Q d 
_ fid + f2 c fj _ 

with Q cS = (fi 2 + 8ft 2 ) 1 / 2 . 

If the solution remains in this state for the whole time, the adiabatic passage from 
initially occupied state a to the target state g takes place provided the Gaussian pulses 
fld(t), fip(t) arrive in a counter-intuitive sequence. 

We are now interested in the linear stability of the dark state (83). To this end, we 
suppose that the solution of (81) evolves in the close neighbourhood of the dark state 

i-e. we express it as a sum 

*(*) = *o(*) + (84) 
where 5&(t) is a deviation of the current solution from the dark state. Inserting this 
expression in (81) yields 

i*o(0 + = /(*o(f) + (85) 

Using (82), one finds 

i^-SV = M5^ - i* - (86) 

where M is a matrix with the elements = J^- with i,j = a, e, g. Here the partial 
derivatives are calculated at the dark state. In system (86) the non-linear terms have 
been omitted. 



1/2 2fi 



^ = < = (83) 
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Specifically, for the system (77,78,79), the matrix M reads 
/ fip^f \ 

m= n P 4>l -(A + i 7 ) n d /2 . (87) 

^ fid/2 i 

(See also (7) in [23]). It is similar to the Hamiltonian (9) of the linear lambda system. 
The difference is that the elements M ae and M ea contain the component of the dark 
state ip® due to the nonlinearity. 

Denoting the Jacobian as A = — iM, we can rewrite the linearised equation (86) as 

—5* = A5V - * - (88) 

In this system, the real parts of the eigenvalues of the Jacobian A determine the stability 
of the dark state. If the matrix M has eigenvalues ui, the Jacobian A has the eigenvalues 
A = — iui. In analogy to the linear lambda system, the eigenvalues oj can be found from 
the characteristic equation 

det||M-/w|| = 0. (89) 
Solving (89) with (87), we find that one root is always zero: 

wi = 0. (90) 
The other two eigenvalues obey the quadratic equation: 

co 2 + (A + i 7 ) w - (S% + 4^ a °| 2 )/4 = 0. (91) 
The corresponding eigenvalues u 2 and w 3 obey the following 

UJ2 + UJ3 = ~ (A + i 7 ), (92) 

oj 2 u 3 = -(^ + 4fi^ a T)/4. (93) 
By setting A = 0, we get the solutions of quadratic equation: 

^2,3 = [-i7 ± (-7 2 + nS + 4nJ|^l 2 ) 1/a ]/2. (94) 

(See also the equations under (7) in [23]). 

In figure 7 we show the dynamics of Re(A(t)) for the non-linear lambda system. 
One can see that this picture reproduces the same behaviour as the corresponding 
dependence for the linear lambda system shown in figure 2. This means that the real 
parts of eigenvalues of the Jacobian for the nonlinear lambda system are the same as 
those for the linear system in the corresponding time intervals: t < t 1 , t > t 2 and 
h < t < t 2 . 

We now adiabatically eliminate the excited state by setting ip e = 0. From (78) we 

get 

^ e = 2(A + i 7 ) ( " p ^ + ^ s) ' (95) 
Inserting this result into (77,79) we obtain a second order system 

r0 a = fi p ^ e , (96) 
r0 g = ifi d ^ e . (97) 
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with ip c given by (95). 

We now discuss the validity of 3D and 2D systems. The 3D system (77,78,79) 
is valid for all times. In the ranges t < t\ and t > t<i the 2D system (96,97) can be 
applied since there are two zero real parts Re(Ai i 2) = 0, and one negative real part 
Re(As) ~ —7. In the range t\ < t < t 2 the process is ID since there is only one zero 
real part Re(Ai) = and two negative real parts Re(A 2i3 ) = —7/2. However, here we 
do not have any ID equation, one can only propose the 2D system (96,97). The search 
for a ID system is a challenging problem. 
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Figure 7. Dynamics of real parts of eigenvalues of the Jacobian for non-linear lambda 
system, computed by (94), and (90). The dashed vertical lines set the boundaries for 
the ID and 2D processes. Here t\ — 0.76 and ti = 5.71. The parameters are as follows: 
A = 0, 7 = 2.0, O = 300.0, t p = 3.8, and t d = 3.0. 

In the figure 8 we have plotted the relevant dynamics for the case of adiabatically 
reduced non-linear lambda system. In the figure 8(a) we show the Gaussian pulses that 
are ordered counter-intuitively. From figure 8(b) one may conclude that the solutions 
of adiabatically reduced system are in good quantitative agreement with the solutions 
of exact system (we do not distinguish them in the figure). In figure 8(c) we see that 
the difference for P c between exact and approximated solutions is significant. It can 
be explained by the fact that the magnitude of probability P c is small. On the other 
hand, the difference in the case of P a , P g is of the same order but we do not distinguish 
it, since the magnitudes of these quantities are much larger. In figure 8(d) we have 
plotted the dynamics of the difference between the exact (P a ) and adiabatically reduced 
(P a d ) solutions of the population P a . (See also the blue line in figure 2(c) in [25]). The 
difference is of the same order as that in figure 8(c). 

In figure 9 we have plotted the dynamics of differences between the populations of 
the current state and corresponding dark state. The difference for the initial ground 
state a deviates up to 0.03 at the end of the passage. Whereas the difference for the 
excited state remains much less. The corresponding difference for state g is almost 
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symmetric to that of the state a in respect to the zero (not shown here). We conclude 
that the process is adiabatic since the solution remains in a close neighborhood to the 
dark state. 




Figure 8. The dynamics (a) of Gaussian pulses (6,7), (b) of populations P a , P g , (c) 
of population P e in enlarged scale, and (d) the difference of P a between exactly (P a ) 
and adiabatically (P^ d ) obtained solutions. In (b,c) the solid thin lines represent the 
solutions of exact (77,78,79), and the (solid thick/dotted) lines are the solutions of 
adiabatically reduced system (96,97) with excited state given by (95). The parameters 
are the same as in figure 7. The initial conditions are Va(0) = 1, ^g(O) = V'e(O) = 0. 




Figure 9. The dynamics of the differences of populations between the current state 
and the dark state, for 5P a {t) = \ip & (t)\ 2 - \ip° a {t)\ 2 (solid line), and SP c (t) = 2\^ c (t)\ 2 
(dashed line) for the nonlinear lambda system. The dynamics of ip a (t),tpe(t) arc 
computed from (77,78,79), and V2(*) i s gi ven by (83). The parameters are the same 
as in figure 7. The initial conditions are V'a(O) = 1, V'g(O) = V'e(O) = 0. 
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5. The non-linear tripod 

We consider the atom-molecule transition in ultracold quantum gases via photoassocia- 
tion. This problem has potential applications during the creation of ultracold molecules 
and quantum superchemistry. The underlying physics is closely related to the STIRAP 
and has been widely studied in the context of atomic physics and quantum optics [3, 52]. 

The level structure of the atom- molecule tripod system is shown in figure 3. In 
ultracold atomic systems, the level a denotes the atomic Bose-Einstein condensation 
(BEC), which couples an excited state of a diatomic molecular BEC via the pump field 
fip(t). Such an excited state is represented by high-lying vibration levels of the single 
excited molecule. The excited state is coupled to the two ground states of the molecular 
BEC, gl and g2, with the strengths fi dl (t) and f2 d2 (t), respectively. 

Assuming a two-photon resonance condition, the four-level non-linear Hamiltonian 
for non-linear tripod takes the form: 

h = -h(A + h)^J c + l(n p ^JJ a + n dl ftj gl + /gg\ 

Here ^> a ,^>t (a = a, e,gl,g2) are the bosonic annihilation and creation operators for 
state a, respectively. To explore the behaviour of the system under time evolution, we 
consider the problem under mean-field approximation, which is reasonable for bosonic 
systems when the number of particles is large compared with unity [29, 37]. In this limit, 
the bosonic operators are replaced by c numbers, and the Heisenberg equation leads to 
the following equations of motion for the probability amplitudes: 

r0 a =O p ^i, (99) 

k/> e = - (A + i7)^e + ^^P^a + T^ 1 ^ 1 + 

+ ^JWg2, (ioo) 

i^ g l = i^dl^e, (101) 

i^ g2 = ifi d2 ^ e . (102) 

The nonlinear term enters here when the molecules are obtained via associating cold 
atoms. 

The normalization reads: 

|^ a (t)| 2 + 2[|^ g i(t)| 2 + |V g2 (*)| 2 + \Mt)\ 2 ] < 1, (103) 

where the equality holds for the initial time. 

Like in the case of linear tripod, the Gaussian pulses are given by (40,41,42). 

Similar to the linear case, we define the following state vector of the system: 
* = [V'ai V'ej V'gi) i J g2] T - The dynamic equations (99,100,101,102) can be rewritten in 
a vector form given by (81), where / is now the vector of non-linear functions on the 
r.h.s of (99,100,101,102). The manifold of the steady states of this system represents the 
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M = 



(107) 



dark state * = [V'ai ^gi> ^g2] T - This manifold has to satisfy (82) and the condition 
of normalization. Since ip® = 0, the dark state obeys the following: 

Q p (^l) 2 + n dl if>° gl + fi d2 ^ g ° 2 = 0, (104) 

i^aT+2i^ g o 1 r+2i^ g o 2 r = i, (los) 

ij° e = 0. (106) 

We are interested in the linear stability of this dark state. Therefore we suppose that 
the solution of (99,100,101,102) evolves in the close neighbourhood of the dark state 
i.e. we express it as a sum given by (84) where S^(t) is the deviation of the current 
solution from the dark state. Hence one arrives at a linearised equation similar to the 
one for the non-linear lambda system (86) where the matrix M reads: 

I ft P v>°* o \ 

IVAa -[A + i 7 ] ^d 2 /2 

o n dl /2 o o 
\ o n d2 /2 o o j 

Note that this matrix is very similar to the Hamiltonian (43) for the linear tripod. 
The main difference between them is dependence of M on ijj®, that arises due to the 
nonlinearity. On the other hand, this matrix is also similar with corresponding matrix 
for the non-linear lambda system. The Jacobian of the linearised system is given by 
A = —\M. The eigenvalues u of the matrix M correspond to eigenvalues A = — iu of 
Jacobian. Solving the eigenvalues problem for matrix M, we get two zero eigenvalues: 

wi, 2 = 0. (108) 
The other two eigenvalues can be found from the quadratic equation 

u? + (A + i 7 ) W - + Q 2 d2 + 4^ a °| 2 )/4 = 0. (109) 
The eigenvalues ou^ satisfy the condition 

+ = - (A + i 7 ), (110) 

u 3 u 4 = - (Q 2 dl + n 2 d2 + 4Q^T)/4- (HI) 
We again assume that A = 0, thus obtaining the following solutions: 

0,3,4 = [-i 7 ± (-7 2 + n2i + ^d 2 + 4^|^ a | 2 ) 1/2 ]/2. (112) 

In figure 10(b) we plot the dynamics of eigenvalues of Jacobian for the non-linear tripod. 
The way of finding the eigenvalues is discussed below. As we saw above, the behaviours 
of corresponding eigenvalues for the linear and non-linear lambda systems was the same. 
Comparing the figures 10(b) and 4 we see that here one can make an identical conclusion: 
the roots behave in the same manner for the linear and non-linear tripods. 

Now we discuss the computing the dynamics of the real parts of eigenvalues, i.e. 
Re(A) = Re(A(t)). The matrix M in the present case depends on ip®, (see M ae and M ea in 
(107)). In the case of non-linear lambda system, the dark state was uniquely defined as 
a function of Rabi frequencies, (83). However, in the present case, for non-linear tripod, 
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the dark state is a manifold that is given by (104,105,106). But we need a definite 
function of time vp® = ijj®(t) in order to get the dynamics of eigenvalues. Therefore we 
use the parametrization of the dark state that was derived in [32]. If the solution of 
(99,100,101,102) evolves on the dark state manifold, we may express that solution in 
terms of only two variables (parameters), [ui(t), U2(t)]: 

1 1/2 



€ = 



4 



u 



1/2 



(113) 



_cos(0)_ 

Mi sin (0) -m 2 cos(0), (114) 
— Mi cos (0) — u 2 sin (0). (115) 

(See the system of equations before (17) in [32]). Here S p = Qdi/Qp and is defined by 
tan (0) = Q^/^di- I n [32] it was also shown that in the adiabatic limit, the parameters 
should obey the equations (see (17) in [32]) 



u l + Q U2 = 0, (116) 

^ ( 1 + A ^ 77^ ) ~ ©«i + 4 
V 4w 2 cos (0) ) dt 

We integrate the system (116,117) and insert its solution in the parametrization 
(113) thus obtaining the necessary dynamics of ip®(t). After inserting this dynamics 
in (112) we get the dynamics of eigenvalues of the Jacobian for the non-linear tripod 
system. 

We may also suppose that the deviation for the amplitude ip a is almost zero, 
8ipa,(t) ~ 0. We thus can make a substitution in (110,111) and (112): 

V>a^a- (H8) 

Subsequently one can numerically solve the system (99,100,101,102). By inserting 
ip a (t) in (110,111) and (112), one gets the approximate dynamics of the eigenvalues. 
Actually, this approach means the analysis of the stability of the current solution 
*(t) = [^ a (t),^ e (t),^ g i(t),^ g2 (t)] T . The dynamics of \ip e (t)\ and |^ c (t)| 2 is plotted 
in figure 10(a). The first dynamics indicates that the magnitude of the r.h.s. of (100) is 
of the order of 0.06. The second dynamics shows that the excited level remains almost 
unpopulated throughout the passage. We therefore conclude that the process is almost 
adiabatic, and one may justify the substitution (118). 

In figure 10(b) we plot the dynamics of Re(A(t)) computed by the both ways. 
The solid line shows the dynamics of Re(A(t)) computed by using the exact value of 
V'a(O) an d l ne dotted line displays the approximate dynamics that is obtained by using 
the substitution (118). We can see from figure 10(b) that the stability of the current 
solution (dotted line) is identical to that of the dark state at the beginning and the 
middle of the process. However, the splitting of the real parts for the approximate 
eigenvalues is slightly delayed with respect to the exact ones. The good quantitative 
agreement of the both results confirms the validity of the approximation (118); it also 
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= 0. 



;ii7) 



Stability of lambda and tripod systems 



25 



shows that the current solution evolves in the close neighbourhood of the moving dark 
state (113,114,115). 
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Figure 10. (a) Dynamics of |^> e (*)| (dotted), and |V> e (i)| 2 (solid) for non-linear tripod, 
both computed by (99,100,101,102); (b) dynamics of real parts of eigenvalues of the 
Jacobian for non-linear tripod, computed by (112), and (108). The solid lines are the 
exact roots, and the dotted lines correspond to stability of the current solution (see in 
the text for details). The dashed black vertical lines set the boundaries for the 3D and 
2D processes. Here t\ — 6.26 and t 2 = 12.29. The parameters are as follows: A = 0, 
7 = 2.0, il = 60.0, t p = 10.7, i d i = 10.0, t A2 = 8.5, K x = 0.75, and K 2 = 5.0. The 
initial conditions are Va(0) = 1, V'gi(O) = ^2(0) = V'e(O) = 0. 



Exactly as in the previous sections, we adiabatically eliminate the excited state by 
setting ip e = 0. From (100) we get 

1 



-(fip^ a + Odi^gl + ^d2^ g 2)- 



2(A + i 7 ) 

Inserting this expression in (99,101,102), we obtain 

k/> a = fip^*^ e , 
i^gl = -fidl^e, 



2 
1 



g2 



(119) 

(120) 
(121) 

(122) 



In these equations we use the expression of ijj e given by (119). 

In figure 11 we have depicted the dynamics of the deviations of the current 
populations from those of the dark state. One can see that all the three populations 
deviate up to 0.05 showing that the solution remains in a proximity to the dark state 
manifold. Therefore one may conclude that the process is adiabatic. 

In figure 12 we have plotted the dynamics for the case of non-linear tripod. Figure 
12(a) shows a sequence of Gaussian pulses. In figure 12(b) we show the dynamics of 
populations of the levels. The solutions of approximated system are in good quantitative 
agreement with those of the exact system (in the figure we do not distinguish them). 
From figure 12(c) we see that the population of the excited state is reproduced by the 
approximated system with a significant error. Again, as in the previous section, we 
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explain this fact by small magnitude of quantity P c . In figure 12(d) we have plotted the 
difference between exact (P a ) and adiabatically approximated {P^ d ) population P a . It 
is of the same order as in the case of excited state (see figure 12(c)). 



6. Some remarks about the one-photon detuning 

Up to now we have been setting A = for all considered systems. In this section 
we shall explore a behaviour of the system in the presence of the non-zero one-photon 
detuning A^O. In that case, when solving the quadratic equations for the eigenvalues 
of the Jacobians, one gets the complex- valued discriminants, D = | | e 1 ^, where \D\ is 
their real amplitude, and f is their phase. The eigenvalues of the Jacobian then read: 



^3,4 — 



+i 



■- 7 ±Vr /2 sm 
2 ' 2 1 1 



+ 



2 2 1 1 



(123) 



Here the real and imaginary parts of the discriminant are given by 

Re(D) = A 2 - 7 2 + fidi + fi d2 + 4^ a T, (124) 
lm(D) = 2 7 A. (125) 

These equations are valid for the non-linear tripod; for the other systems we get the 
similar expressions. If A = (as assumed previously), the imaginary part becomes zero, 




Figure 11. The dynamics of the differences of populations between the current state 
and the dark state for the nonlinear tripod, SP a (t) = \ip a (t)\ 2 — IV'a (*)| 2 (dashed line), 
SP c (t) = 2|VcW| 2 (dotted line), SP gl (t) = 2|V> g i(t)| 2 - 2\iP° gl (t)\ 2 (solid line), and 
8P&(t) = 2\Tp g2 (t)\ 2 -2\ip° 2 {t)\ 2 (dash dotted line). The dynamics of tp a (t), tp e (t), ip s i(t) 
and ip g 2(t) are computed from (99,100,101,102), and V'a(*)>^'gi(*)>^'g2(*) are found from 
(113,114,115) and (116,117). The parameters are the same as in figure 10. The initial 
conditions are Va(0) = 1, Vgi(0) = ^g2(0) = V'e(O) = 0. 



Stability of lambda and tripod systems 



27 



300- 



- 200- 



— 



Q_ 



(a) 


5 10 


15 


20 






t 






1.5 










V 

1 


1.0 






i 






0.5- 










o 


0.0- 








A 




(c) 


) 


5 10 
t 


15 


20 




Figure 12. The dynamics (a) of Gaussian pulses (40,41,42), (b) of populations P a , 
Pgi, Pg2, and (c) of population P c in enlarged scale. In (d) it is shown the difference 
between the exact (P a from (99,100,101,102)) and adiabatically reduced (P a d from 
(120,121,122)) populations in enlarged scale. In (b,c) the solid thin lines represent the 
solutions of exact (99,100,101,102), and the (solid thick/dotted) lines are the solutions 
of adiabatically reduced system (120,121,122) with excited state given by (119). The 
parameters are as follows: A — 0, 7 = 2.0, fi = 60.0, t p = 11.5, idi = 10.0, id2 = 8.5, 
Ki = 0.75, and K 2 = 5.0. The initial conditions are the same as in figure 10. 



i.e. lm(D) = 0; the phase may be either or 7r; in the interval t± < t < £2, it is (p = 0, 
and in the ranges t < ti, t > t<i it is ip — ir. For ip = we have Re(A 3i4 ) = —7/2, 
and for <p = ir we get Re(As) ~ 0, Re(A4) ~ —7. We have got these results for all the 
considered systems (see e.g. figure 10(b)). However, in the case of the non-zero one- 
photon detunings, in the interval t x < t < t 2 , these real parts are no longer coinciding; 
they are symmetrically surrounding the value —7/2, and the difference between them 
becomes equal to 

ReA 3 - ReA 4 = |/J| 1/2 sin (|) ~ 

~ Im ^) (126) 
" 2[Re(£>)]V2- 

The latter approximation is valid for the small values of lm(D) << R,e(D). Such a 
situation takes place in the middle of the passage, when the Rabi frequencies are large 
compared to the one-photon detuning and losses. 

We thus conclude that for such small detunings the difference between the negative 
real parts remains small, and our statements about the reduction of dimension remain 
valid. 
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7. Conclusions 

We have analysed the adiabatic reduction of the dimension of the linear and non-linear 
three- and four-level systems. By evaluating the corresponding Jacobians and computing 
the dynamics of real parts of their eigenvalues (the non-zero eigenvalues are found from 
quadratic characteristic equations), one may define the dimensionality of the processes. 
This dimensionality is given by the number of zero real parts since the negative real 
parts cause the contraction of the nearby solutions towards the dark state. At the 
beginning and the end of the dynamics, there is always only one negative real part. 
Hence one may eliminate only one state representing the excited state. In the middle 
of the process, one of the zero real parts becomes negative thus making the number of 
negative real parts equal to two. In this time interval we may eliminate two variables 
corresponding to the excited and bright states respectively. For linear systems, we have 
eliminated both excited and bright states. However, for non-linear systems, we have 
restricted ourselves by eliminating the excited state. This is due to the fact that the 
definition of a bright state for the non-linear systems is not available. We suppose that 
the remaining stable degrees of freedom in the non-linear systems can be eliminated by 
using the asymptotic methods of non-linear dynamics. 

The main finding of this work is revealing that the whole STIRAP evolution for 
all considered systems is divided into three time intervals with different number of the 
negative real parts of the Jacobians. The evolution of the real parts is equivalent for the 
corresponding linear and non-linear systems (as one can see in figures 2,4,7, and 10(b)). 
This suggests that the non-linear systems may be potentially reduced as the linear ones. 
Physically this means that the considered three/four- level schemes may be regarded as 
schemes with lower dimension, i.e. with fewer levels involved. In the time intervals 
t < ti and t > t 2 , the initially three- level system is effectively a two- level one, and in 
the range of t\ < t < t 2 it contains a single level. Analogously, the initially four-level 
scheme may be regarded as a three-level (two-level) system in the time ranges, where 
t < t±, t > t 2 (h < t < t 2 ). 

A sensitive problem is the definition of a dark state for the non-linear tripod. In 
the case of non-linear lambda system, the dark state is a moving point (83) in the phase 
space. However, for the non-linear tripod we have a manifold (104,105,106) of dark 
states. If one wishes to get the dynamics of real parts of eigenvalues of the Jacobian, 
one needs a definite value of complex amplitude i/j® belonging to the manifold. We here 
use two ways for the stability analysis of the dark state. The first way is to parametrize 
the dark state manifold by using the method developed in [32]. This method enables 
one to find the definite dynamics of ip®{t). We thus managed to find the exact dynamics 
of eigenvalues. The second way is to simply substitute using (118) the value h Y 
the current solution ip a (t) that is found from underlying equations (99,100,101,102). 
Actually, the substitution (118) means we are investigating the stability of the current 
solution instead of that for the dark state. In fact, figure 10(b) shows that the real parts 
of the eigenvalues evolve almost identically. The only difference is that the splitting of 
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real parts for approximate eigenvalues is slightly delayed. Such a coincidence shows that 
the current solution evolves in the close neighbourhood to the motion of the parametrized 
dark state (113,114,115). It is also to be noted that the magnitude of |^ e | is always small, 
and the excited state remains almost unpopulated as we can see in figure 10(a). 

It is noteworthy that a related approach was used in [30], where a feedback 
control scheme was presented that designs time-dependent laser-detuning frequency 
to suppress possible dynamical instability in coupled free-quasibound-bound atom- 
molecule condensate systems. It was proposed to perform a substitution analogous 
to (118) which was used for solving the control problem. On the other hand in our work 
this substitution was made for the stability analysis of the dark state. 

It is also important to note that in the lambda and tripod systems, we have 
phenomeno logically included the loss coefficient 7. This was done by making the one- 
photon detuning to be a complex number, i.e. by replacing A — > A + 17. Here A 
is again a one-photon detuning, and 7 determines the losses. In our work, we have 
considered the cases where A = and 7 > 0, i.e. the one-photon resonances with 
losses. We stress that the presence of non-zero losses 7 makes the adiabatic reduction 
easier to implement. The losses cause the appearance of two negative real parts of 
eigenvalues of the corresponding Jacobians. On the other hand, it was shown that the 
losses decrease the transfer efficiency [8] , which decreases exponentially with the (small) 
decay rate. However the range of decay rates, over which the transfer efficiency remains 
high, appears to be proportional to the squared pulse area. Hence, by choosing high 
pulse areas one may preserve the high transfer efficiency 

Another question is a possible presence of the one-photon detuning in the considered 
processes. As it was shown in section 6, the relatively small one-photon detuning does 
not alter our conclusions about the reduction of dimension in the three- and four-level 
systems considered here. This happens if the Rabi frequencies are large compared to 
the one-photon detuning and loss rates. 
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